1. Field of the Invention
The present invention relates to modeling of semiconductor devices and specific physical phenomena.
2. Description of the Related Art
The basic equations applied to modeling of semiconductor devices are represented generally in the following form: EQU .differential.p/.differential..tau.=(-1/q)(divJ.sub.p)+G-U (1) EQU .differential.n/.differential..tau.=(1/q)(divJ.sub.n)+G-U (2) EQU J.sub.p =-qD.sub.p (gradp)-q.mu..sub.p p(grad.psi.) (3) EQU J.sub.n =qD.sub.n (gradn)-q.mu..sub.n n(grad.psi.) (4) EQU div(grad.psi.)=(-q/.epsilon.)(N.sub.d -N.sub.a +p-n) (5)
These five equations are applied commonly to the conventional modeling method and the modeling method according to this invention. Equations (1) and (2) are a continuity equation for holes and a continuity equation for electrons, respectively. In these equations, .tau. is physical time, and G and U denote the generation and recombination of carriers. Equations (3) and (4) explicitly describe a hole current density and a electron current density, respectively. As can be understood from equation (3), the hole current density J.sub.p is an algebraic sum of a diffusion current component (the first term) proportional to the carrier density gradient gradp and a drift component (the second term) proportional to the potential gradient grad .psi. and the carrier density. Also, as can be seen from equation (4), the electron current density J.sub.n is an algebraic sum of a diffusion current component (the first term) proportional to the carrier density gradient gradn and a drift component (the second term) proportional to the potential gradient grad .psi. and the carrier density. Equation (5) is a Poisson's equation representing the relationship between the potential .psi. and the space charge density.
The conventional method of modeling semiconductor devices will be first explained. In the first step of obtaining the numerical solutions of equations (1) to (5), these equations are rewritten from continuous system (i.e., differential equations) to discrete system (i.e., difference equations). To rewrite equations (1) to (5) so, it is necessary to divide the device space, which is to be analyzed, into a rectangular meshpoint array. Mesh division for a two-dimensional space model can be schematically illustrated in FIG. 1. Let us define the meshpoint spacing hx in the x-direction and the meshpoint spacing hy in the y-direction as is illustrated in FIG. 1, and also let us define the midpoint spacing h.sub.x ' in the x-direction and the midpoint spacing h.sub.y ' in the y-direction as follows: EQU h.sub.x '(M)=(1/2)[h.sub.x (M'-1)+h.sub.x (M')] (6-1) EQU h.sub.y '(N)=(1/2)[h.sub.y (N'-1)+h.sub.y (N')] (6-2)
From equations (3), (6-1), and (6-2), the hole current divergence div J.sub.p can be represented by the following approximation in the discrete system: ##EQU1## As for the electron current, the divergence div J.sub.n can be rewritten in the same way. Therefore, equations (1) and (2) are rewritten to: ##EQU2##
In equations (8-1) and (8-2), the carrier-generation term G is omitted for simplicity's sake, since this term is not required in analyzing semiconductor devices under ordinary operating conditions. Also for simplicity's sake, the physical time differential terms .differential.p/.differential..tau. and .differential.n/.differential..tau. are set at zero, on the assumption that calculations are performed to find a DC steady-state solution. It should be noted, however, that calculations can be performed essentially in the same way to obtain a non-steady-state solution, where the physical time differential terms are not zero.
Likewise, Poisson's equation (5) can be expressed as follows: ##EQU3## where, for the sake of simplicity, doping function, i.e., the difference .GAMMA. between the donor concentration and the acceptor concentration, is set as .GAMMA.=N.sub.d -N.sub.a.
Here, the carrier recombination term U(M,N) is expressed in the Shockley-Read-Hall form for the ensuing development. Namely: ##EQU4##
Next, equations (3) and (4) for currents are discretized, by applying the Scharfetter-Gummel approximation in order to, as is known in the art, ensure the stability in numerical analysis. As a result of this, the x-component and y-component of the hole current density, and those of the electron current density are given as follows: EQU J.sub.px (M')=[q/h.sub.x (M')][.lambda..sub.px1 (M')p(M N)+.lambda..sub.px2 (M')p(M+1,N)] (11-1) EQU J.sub.py (N')=[q/h.sub.y (N')][.lambda..sub.py1 (N')p(M N) +.lambda..sub.py2 (N')p(M,N+1)] (11-2) EQU J.sub.nx (M')=[q/h.sub.x (M')][.lambda..sub.nx1 (M')n(M N)+.lambda..sub.nx2 (M')n(M+1,N)] (11-3) EQU J.sub.ny (N')=[q/h.sub.y (N')][.lambda..sub.ny1 (N')n(M N)+.lambda..sub.ny2 (N')n(M,N+1)] (11-4)
The .lambda.-items in equations (11-1) to (11-4) are defined by the following equations: EQU .lambda..sub.px1 (M')=[.mu..sub.p (M')/.theta.(M')][.beta.(M')/(1-e.sup.-.beta.(M'))] EQU .lambda..sub.px2 (M')=[.mu..sub.p (M')/.theta.(M')][.beta.(M')/(1-e.sup..beta.(M'))] EQU .lambda..sub.py1 (N')=[.mu..sub.p (N')/.theta.(N')][.beta.(N')/(1-e.sup.-.beta.(N'))] EQU .lambda..sub.py2 (N')=[.mu..sub.p (N')/.theta.(N')][.beta.(N')/(1-e.sup..beta.(N'))] EQU .lambda..sub.nx1 (M')=[.mu..sub.n (M')/.mu..sub.p (M')].lambda..sub.px2 (M' ) EQU .lambda..sub.nx2 (M')=[.mu..sub.n (M')/.mu..sub.p (M')].lambda..sub.px1 (M' ) EQU .lambda..sub.ny1 (N')=[.mu..sub.n (N')/.mu..sub.p (N')].lambda..sub.py2 (N' ) EQU .lambda..sub.ny2 (N')=[.mu..sub.n (N')/.mu..sub.p (N')].lambda..sub.py1 (N')(12)
Thus, basic equations (1) to (5) for modeling semiconductor devices have been rewritten to discrete ones. The equations which we are to solve are: equation (9); equation (8-1) with equations (11-1) to (11-4), auxiliary relations (12) and recombination equation (10) substituted therein; and equation (8-2) with equations (11-1) to (11-4), auxiliary relations (12) and recombination equation (10) substituted therein. These three equations contain three unknown quantities p, n, and .psi.. However, it is difficult to solve equations (8-1) and (8-2) unless these equations are transformed, since these equations each contain nonlinear relations for the three unknown quantities. Hence, equations (8-1) and (8-2) are expanded in Taylor series for the unknown quantities, neglecting quadratic terms and higher-order terms, thereby to derive linear equations for the changes .delta.p, .delta.n, and .delta..psi. of the unknown quantities p, n, and .psi..
More specifically, we first assume that the unknown quantities are each the sum of a trial value (assigned with superscript 0) and a change. That is: EQU p(M,N)=p.sup.0 (M,N)+.delta.p(M,N) EQU n(M,N)=n.sup.0 (M,N)+.delta.n(M,N) EQU .psi.(M,N)=.psi..sup.0 (M,N)+.delta..psi.(M,N) (13)
Then, current equations (11-1) to (11-4) are expanded in Taylor series, neglecting the quadratic terms and the higher-order terms, thus yielding the following results: ##EQU5## where superscript 0 is assigned to quantities which correspond to p.sup.0, n.sup.0, and .psi..sup.0.
In addition to the current equations, recombination term U(M,N) is also nonlinear. Recombination term U(M,N) is therefore expanded in Taylor series. That is: EQU U(M,N)=U.sup.0 (M,N)+U.sub.p.sup.0 (M,N).delta.p(M,N)+U.sub.n.sup.0 (M,N) .delta.n(M,N)
where ##EQU6##
Equations (13), (14-1), (14-2), (14-3), (14-4), and (15) are substituted in equations (8-1), (8-2) and (9), yielding the following matrix-vector equation for the two-dimensional problem: ##EQU7## where, by definition of matrixes and vectors: ##EQU8## The matrix elements and the vector elements can be more explicitly represented as follows: ##EQU9##
As can be understood from the above, matrix-vector equation (16) can be derived by linearizing discrete equations (8-1), (8-2) and (9).
Equation (16) corresponds to one of many meshpoints arranged in the x-direction and the y-direction. Hence, in order to obtain a solution for the whole device plane, equations (16) for the number of meshpoints must be solved simultaneously. In the case where the numbers of meshpoints forming the whole device plane are in the ranges of, for example, 1.ltoreq.M.ltoreq.K and 1.ltoreq.N.ltoreq.L, the whole equations (16) will be represented as FIG. 2. For simplicity's sake, A.sub.T (M,N), .THETA..sub.T (M,N), and F.sub.T (M,N) have been rewritten to A.sub.MN, .THETA..sub.MN, and F.sub.MN, respectively.
In FIG. 2, each matrix, for example B.sub.11, has 3.times.3 dimensions. Therefore, the matrix system defining the whole device plane has (3KL).times.(3KL) dimensions. In the case of a standard device model, L=30 and M=40. Hence, the matrix system for the standard device model has the dimensions of 3,600.times.3,600=12,960,000=1.296.times.10.sup.7.
Various methods are available which can be used to solve the matrix-vector equation schematically represented in FIG. 2. The most commonly employed method is Gaussian elimination. Gaussian elimination method involving a matrix system of the size specified above can be performed by means of a high-speed, large-memory computer now commercially available. Other alternative methods that can be employed to solve the matrix-vector equation of FIG. 2 are iteration methods which are effective particularly for the numerical analysis of very large matrix systems.
A matrix-solving method is employed not only in device modeling, but also in solving physical and engineering problems. Since the method is a calculation process well known in the art, it is assumed here that equation (16) has been solved for all space of the device. Obviously this means that the solutions of the linearized equation (16), which originates from the nonlinear equations (8-1), (8-2) and (9), have been obtained. If these solutions are substituted in the nonlinear equations (8-1), (8-2) and (9), both sides of the equation will not be equal in most cases, due to nonlinear effect. Let us therefore substitute the solution [p(M,N), n(M,N), .psi.(M,N), 1.ltoreq.M.ltoreq.K, 1.ltoreq.N.ltoreq.L] as a set of trial values [p.sup.0 (M,N), n.sup.0 (M,N), .psi..sup.0 (M,N), 1.ltoreq.M.ltoreq.K, 1.ltoreq.N.ltoreq.L], and then execute the same calculation process as mentioned above. Further, let us repeat this calculation process, a necessary number of times, then we will obtain solutions which completely satisfy equations (8-1), (8-2) and (9).
The matrix-solving method described above can be achieved to solve very large problems by means of the computer now available which has a great memory storage and can perform operations at high speed. However, the higher the performance of the computer, the larger a problem people wish to solve by means of computer. In view of this, it is still one of the greatest difficulties to solve matrix problems at a sufficiently high speed in order to accomplish successful modeling of semiconductor devices.
To solve these equations at high speed, methods and apparatus have been proposed which can solve particular equations, without applying matrix-solving methods, thereby to analyze the operations of semiconductor devices and specific physical phenomena. The present invention relates to an improved method of determining the sensitivity coefficients or the amplification gains, which helps to enhance the efficiency of the analysis of the semiconductor device operations and the specific physical phenomena.